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We present a theoretical study of indirect exchange interaction between magnetic adatoms in 
graphene. The coupling between the adatoms to a graphene sheet is described in the framework of 
tunneling Hamiltonian. We account for the possibility of this coupling being of resonant character 
if a bound state of the adatom effectively interacts with the continuum of 2D delocalized states 
in graphene. In this case the indirect exchange between the adatoms mediated by the 2D carriers 
appears to be substantially enhanced compared to the results known from Ruderman-Kittel-Kasuya- 
Yosida (RKKY) theory. Moreover, unlike the results of RKKY calculations in the case of resonant 
exchange the magnetic coupling between the adatoms sitting over different graphene sublattices do 
not cancel each other. Thus, for a random distribution of the magnetic adatoms over graphene 
surface a non-zero magnetic interaction is expected. We also suggest an idea of controlling the 
magnetism by driving the tunnel coupling in and out of resonance by a gate voltage. 

PACS numbers: 68.65.Pq, 75.30.Hx, 75.50.Pp, 75.75.-c 


I. INTRODUCTION 

Magnetism in graphene has been attracting quite a 
lot of interest, probably no less than any other physical 
property of this novel material. In particular, quite a few 
experiments are focused on introducing magnetic prop¬ 
erties by doping graphene with magnetic atomf®. The 
magnetic adatoms deposited onto graphene surface in a 
moderate sheet density can couple to the graphene and 
participate in the indirect exchange interaction mediated 
by the free carriers available in graphene thus provid¬ 
ing an analog of a dilute magnetic semiconductor. The 
indirect exchange is usually treated with RKKY theory 
which considers the spin-spin interaction between impu¬ 
rities mediated by delocalized carriers. As we have shown 
in our previous works RKKY approach needs to be mod¬ 
ified in the case when the impurity is coupled to the free 
electron gas by means of resonant tunnelling. Basically, 
the spatial separation between the magnetic center and 
free carriers leads to decrease of the indirect exchange due 
to weak wavefunctions overlap. However, if the magnetic 
center posess a bound state with an energy level match¬ 
ing the delocalized carriers energy range, the indirect ex¬ 
change is dramatically increased due to resonant tunnel¬ 
ing effedP'^i This phenomenon was observed in GaAs 
based heterostructures having an In^Gai-^jAs quantum 
well (QW) and Mn <5-layer in the vicinity of Qv#tni. 

In this paper we analyse pair indirect exchange inter¬ 
action between magnetic adatoms via delocalized carriers 
in graphene. For the non-resonant case when the bound 
state energy does not match the 2D continuum energy 
range this approach resides to conventional RKKY the¬ 
ory in graphene with some minor model-specific adjust¬ 
ments. The RKKY interaction in graphene has been in¬ 
tensively studied theoretically 12 It was shown that 
the indirect exchange crucially depends on the magnetic 
centers position with regard to the graphene sublattices 


A and B. If the magnetic centers are located at the 
vertices of the same sublattice the interaction is ferro¬ 
magnetic while for the opposite sublattices or plaquette 
configuration it is antiferromagnetic. For the undoped 
graphene the interaction energy decreases aJHHHE7] i j p>'i 
with distance R between the impurities, while for doped 
graphene the dependence i d — ■*— I 1/R 2 . For the case of 
magnetic adatoms being weakly coupled to the graphene 
the RKKY interaction is simply damped by a factor 
of exp (—4d/do) with d being the distance between an 
adatom and the graphene and do characterizing decay of 
the graphene wavefunction in the direction normal to its 
surface. On the contrary, if an impurity has a resonant 
bound state the situation changes dramatically. In the 
resonant case the energy of the indirect exchange inter¬ 
action is strongly enhanced and the interaction type (fer¬ 
romagnetic or antiferromagnetic) can be the same both 
for AA and AB adatom configurations as will be shown 
below. 


II. THEORY 

Let us consider two magnetic adatoms located on a 
graphene sheet as shown in Fig.l. Magnetic adatoms 
have spins labelled as R , I 2 and the bound states with 
energy levels ei, 62 - We assume a tunnel coupling which 
leads to a hybridization of the adatoms bound states de- 
screte energy levels and graphene continous spectrum. 
Such hybridization is described by the well-known Fano- 
Anderson modePP For a non-resonant case when the the 
discrete level energy lies outside of the occupied contin¬ 
uum spectrum the perturbation theory of indirect ex¬ 
change (RKKY) can be applied using the hybridized 
wavefunctions as the ground state. However, if the dis¬ 
crete level energy lies within the 2D continuum the inter¬ 
action becomes resonant, i.e. the substantial part of the 
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hybridized wavefunction resides at the magnetic centers. 
At that, even a small perturbation of the bound state 
reflects in a strong change of the hybridized wavefunc¬ 
tion. This effect does not allow one to apply the RKKY 
theory which technically leads to a divergence of the en¬ 
ergy correction due to resonant denominators in the sec¬ 
ond order perturbation expressions. An alternative ap¬ 
proach, which avoids using the perturbation theory was 
suggested in our previous works- 7 . The approach still 
assumes weak exchange coupling regime so the magnetic 
adatoms spins are treated as classical magnetic moments. 



FIG. 1: (Color online) The illustration of indirect exchange 
interaction of adatoms located on graphene. Two adatoms 
have spins I \, I 2 and may be coupled to the bought graphene 
sublattice. 


positions of the two adatoms. We express the whole 
Hamiltonian (JT]) in the second quantization as follows: 

H 0 =£ifi fi + £ 2 fth + v F 

P 

- \p\dp K d P K + |p| CpK'CpK' - \p\d+ K ,d P K>), 

H t =t 1A d^Ji +t 1B b^Ji +t 2 Aa B2 h +t2B~b B2 f 2 +h.c ., 
Hj=X 1 f+f 1 +X 2 f/f 2 . (4) 

Here £i >2 are localized states energies, 
Cp k i dpK > bpK 1 : dpK' > bpK 1 dpK , CpK', d p x r are the 
creation and annihilation operators of electrons (c) and 
holes (d) with discrete momentum p at K and K'. 
points respectively. Vf is the Fermi velocity in graphene, 
Ai j2 = A,S7i, 2 |(/?(0)| 2 , where S, I\ t2 are the electron and 
the adatoms spin projections respectively, <^(0) is the 
bound state wave function at the magnetic adatom site. 

To account for the coupling of each adatom to 
the graphene we introduced two complex parameters 
i = 1,2, describing the tunnel coupling of ith 
adatom to the A and B graphene sublattice, respec¬ 
tively. With that the tunneling part Ht hr W) is writ¬ 
ten in real space representation, a/./b/. being the cre¬ 
ation operators for real-space states at A and B sublat¬ 
tice. Along with the part of the Hamiltonian Hq (( 4 ]) de¬ 
scribing graphene the tunneling part must be rewtitten 
in the momentum representation. Using tight-binding 
graphene Hamiltonian 20 1 operators in real-space are ex¬ 
pressed through those in reciprocal k-space as follows: 


We consider the following Hamiltonian: 


H = Ho + H t + Hj , (1) 

where Ho describes the non-interacting adatoms and 
graphene, Ht is the tunneling term, which is spin inde¬ 
pendent, Hj describes the exchange interaction between 
the adatoms core spins and the graphene electrons which 
have tunneled to the adatom. We assume there is no di¬ 
rect exchange interaction between the adatoms. For the 
weak coupling the adatoms spins act as parameters defin¬ 
ing the potential energy profile, thus the Hamiltonian ([T]) 
becomes different for parallel and antiparallel spin con¬ 
figuration. The difference in the system energy for par¬ 
allel and antiparallel spin configurations of the magnetic 
adatoms can be inter preted as the indirect exchange in¬ 
teraction energy^EESnSl 

E ex = Ef | — .Bp(2) 

We assume the exchange hamiltonian Hj as a contact 
interaction at the adatoms sites: 

Hj = AA(r — Rr)SIi + A£(r — R 2 )SI 2 , (3) 

where S is the electron spin operator, I 12 are the 
adatoms spins, A is the exchange constant, Rp 2 denote 
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( 5 ) 


where R„ is a position of a carbon atom in the sublattice 
A, Op is a polar angle of the momentum vector p (here the 
coordinate system is the same as iiP). We assume that 
the first adatom is located at the the origin Ri = (0, 0) 
and the second adatom is located at R 2 = (x, y) = (r, 0 T ) 
in Cartesian and polar coordinate system, respectively. 

We further proceed with diagonalization of ([l]) on a 
discrete basis in the same manner as described irP. It 
appears to be more convenient to use a cylindrical dis¬ 
crete basis as the boundary conditions for the graphene 
in these coordiantes are simpler to define. Let us put the 
system in a finite cylindrical box with radius L. Then one 
could specify a basis of the pure graphene eigenfunctions 
using any type of boundary conditions. The graphene 
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eigenfunctions in cylindrical basis are: 
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where up and down sign correspond to the electrons and 
holes, respectively. We take the boundary conditions for 
K and K' valleys in a general form as discussed irPH 
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here parameter /3 describes the boundary. Using (J7| we 
get discrete energy levels for I\ and K' valleys: 
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where ” +” corresponds to the holes and ” —” corresponds 
to the electrons, to, n are integer values, to is the quan¬ 
tized momentum amplitude n is the cylindrical harmonic 
number. Now let us consider the total Hamiltonian ([lj, 
we write its eigenfunctions | \H) as an expansion: 


W — ( a l/l + + a 2ft + ^pn^pnK + K™ ^pnK~ 
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where the coefficients cii,b ^ are to be determined. The 
indexes p , n characterize the cylindrical basis states ©• 
The tunnelling matrix elements in the cylindrical basis 
(6) could be obtained from those in the plane waves basis 
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Finally, in the basis 

(fi 10), / 2 + 10), {c+ *10)}, { d+ nK 10)}, {c+ x , 10)}, { d+ nK , 10)}) 
the Hamiltonian 0 reads: 
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The eigenvalues of H (11) can be found analytically, the 
details are presented in Appendix A. The indirect ex¬ 
change energy is the difference between parallel and an¬ 
tiparallel adatoms spin configurations, in order to find 


the total energy difference one should sum up the differ¬ 
ence over all occupied states: 

E ex = E [ E m4n,s)-E m!i (it,s)], (12) 

m,i s=± 1/2 


where to = 0 ,1, 2,... and * = 1,2 characterise the modi¬ 


fied energy levels (eigenvalues of (111) and the outer sum 
is over all occupied energy levels (see Appendix A for de¬ 
tails). Also for each of the adatoms configurations (ft 
denotes parallel spin configuration, fj- is f° r antiparallel 
spin configuration) the summation is done over the elec¬ 
tron spin projection s. The adatoms spin configurations 
and electron spin prohection enter (11) only through 
the parameters Ai j2 , which are of the same magnitude 
I Ai| = |A21 = A but can have different sign, for example 
for ft) s = —1/2 these parameters are Ai = A 2 = —A. 
In the following we keep the leading order in the tun¬ 
nelling parameters. As shown in Appendix B the indirect 
exchange energy contains a fast oscillating factors with 
a characteristic spatial period 2tt/K, that is the length 
of the carbon-carbon bond, the same applies to a con¬ 
ventional RKKY interaction in grapheneASUT Here we 
average over these short-wavelength oscillations and also 
average on the polar angle of describing the adatoms po¬ 
sition (for details see Appendix B). Finally, we get for 
the resonant indirect exchange interaction in graphene: 
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t AA ~ ^4 (tlAt* 2A + tiBt*2 B )(dlAt2A + Vb^b), 

tab = ^4 (Vi^.| 2 |*2b| 2 + |Ub| 2 |*2a| 2 ), 

where t ~ 2.8 eV - is the graphene’s nearest-neighbor 
hopping energj®°( Jj, IV, - are Bessel and Neumann func¬ 
tions of *th order, respectively. 


III. DISCUSSION 


Let us assume for simplicity that the adatoms bound 
states have the same localization energy, £\ = e 2 = £o- 
For the resonant case, when the bound state energy lies 
within the 2D continuum spectra, £q € (—00 ,Ep), from 
(13) the exchange energy can be roughly estimated as: 


Eres ~ — V^ £ 0j T — ( t AA + t Ab) ^ (15) 


To consider the non-resonant case let us assume the local¬ 
ized state energy lying above the Fermi level, £q Ep- 
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In this case the arctangent argument has no poles 
and arctan is to be replaced by its argument (due to the 
smallness of the tunneling parameter). Further calcula¬ 
tions are similar to those discussed iiffl We get: 


E, 


nr ~ 4 

net 


J dE E 2 f(E,r) = 

— OO 

(^f(£) [ “ t1a + 3t1s] ' (16) 


The expression (16) is consistent with the previ- 



FIG. 2: (Color online) Indirect exchange interaction energy 
vs distance between adatoms. Localized state’s energy is eo = 
— 100 meV, Fermi level is Ef = 0 meV. Red curve - adatoms 
located on same sublattice (AA), blue dash curve - adatoms 
located on different sublattice (AB). 


ously k nown resu lts of conventional RKKY theory in 
graphen eP^ -l— with the characteristic dependence 
E nr ~ A 2 T 4 . While in metals the the RKKY-type 
indirect exchange interaction is usually ferromagnetic 
at small distances between the magnetic impurities, in 
graphene it depends on the adatoms position. As follows 
from (16) if both adatoms are located at the same sublat¬ 
tice (AWor BB) the interaction is ferromagnetic while for 
different subllatices configuration (AB, BA) it is antifer¬ 
romagnetic. This is in exact agreement with the RKKY 
theory in graphene 1 ^ However, for the resonant case a 
picture of the indirect exchange interaction appears to be 
radically different from what is expected from the con¬ 
ventional RKKY theory. In order to apply the resonant 
exchange theory to any real case one should be aware of 
typical values of the parameters ente ring (13). Thos e can 
be estimated as A ~ 1 4-10 meV® 1 3, T ~ 0.1 4- 1 eV 1 
The STM experiment^ have shown that Co adatoms 
have localized energy level at £q ~ —100 meV and the 
Fermi level in graphene can be effectively controlled by 
gate voltage in the range Ep ~ —300 4- 300 meV'H 


Figs® show the resonant indirect exchange interac¬ 
tion as a function of the distance between the adatoms 
calculated according to the formula (131 for different 
adatoms configurations (AA and AB). For this calcu¬ 
lation we took the following values: A = 1 meV. For 
the AA configuration the tunneling parameters describe 
the coupling of the adatoms to only one sublattice (A) 
tiA = hA = 1 eV, tip = t 2 B = 0. For the (AB) config¬ 
uration we took tiA = t 2 B = 1 eV, t±B = <24 = 0. We 
analyse two different situations regarding the position of 
the bound state energy related to the Dirac point: (i) 
£q = —100 meV, Ep = 0 meV, (ii) £q = 100 meV, 
Ep = 110 meV. In both cases the resonance condition 
is satisfied. Case (i) is presented in Fig|2] Note that 
in contrast to the conventional (non-resonant) RKKY 
theory the interaction between the adatoms is antifer¬ 
romagnetic at small distances for both AA and AB con¬ 
figurations. The reason is that unlike in the RKKY case 
in the resonant case the tunnelling is enhanced for the 
graphene electrons having the same energy as the local¬ 
ized state and for small distance both terms in (14) ap¬ 
pear to be of the same sign. In the case (ii) the bound 
state energy level lies in the conduction band above the 
Dirac point. The calculated indirect exchange energy is 
presented in Figj3] In this case for both AA and AB 
configurations the interaction is ferromagnetic at a small 
distance. Note, that because the contributions from AA 
and AB configurations do not compensate each other, 
the plaquette configuration or random distributions of 
the magnetic adatoms would also result in a ferromag¬ 
netic interaction. Thus, the resonant exchange can be 
revealed in an experiment with a random distribution of 
the adatoms. For both resonant cases (i) and (ii) dis¬ 
cussed above the characteristic wavelength for the oscil- 



FIG. 3: (Color online) Indirect exchange interaction energy 
vs distance between adatoms. Localized state’s energy is £o = 
100 meV, Fermi level is Ef = 110 meV. Red curve - adatoms 
located on same sublattice (AA), blue dash curve - adatoms 
located on different sublattice (AB). 
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lations is governed by the resonant energy: A ~ a ■ t/e^ 
while in conventional RKKY theory it is related to the 
Fermi energy A ~ Vp/sp- 

An important feature of the resonant indirect exchange 
is that its strength can be several orders of magnitude 
higher compared to the non-resonant case. This effect 
has been described for semiconductor heterostructure^ 
and carbon nanotube^. The physical reason for such 
an enhancement is simple, the resonant tunneling makes 
the free carriers wavefunction effectively penetrates onto 
the magnetic center where it participates in the direct 
exchange. The enhancement although not that large 
has been seen experimentally for hybrid semiconductor- 
Mn heterostructurc 11 . Here we suggest a natural idea 
of modulating the magnetic properties of the system by 
adjusting the Fermi level in order to drive the system 
in and out of the resonant condition. To illustrate the 
effect we plotted the indirect exchange energy vs Fermi 
level position in Fig|4] The distance between the adatoms 
located on same sublattice (AA configuration) was kept 
R = |Ri — R.21 = 2 nm and the bound state energy 
was taken £ 0 = —100 meV. As seen from the Fig. [4] 
the indirect exchange energy changes by several orders 
of magnitude while changing the Fermi level position Ep 
by ~ 100 meV. Such a modulation of the Fermi level po¬ 
sition can be rather easily achieved by applying the ap¬ 
propriate gate voltage. Thus, along with the modulation 
of the interaction range 12 one could effectively control its 
strength using the gate voltage. 



E p (meV) 


FIG. 4: (Color online) Indirect exchange interaction energy 
as a function of the Fermi level position. Bound state energy 
is £o = —100 meV, the distance between adatoms is R = 2 
nm, AA configuration. 


IV. SUMMARY 

In conclusion, we have developed an indirect exchange 
interaction theory for magnetic adatoms on graphene. 


The presence of resonant localized states at the adatoms 
strongly affects the picture of indirect interaction. The 
indirect exchange is strongly enhanced whenever the 
adatom bound state energy level lies within the energy 
range of occupied states in graphene. Thus, the indi¬ 
rect interaction strength can be effectively controlled by 
adjusting the Fermi level position in graphene with an 
external gate voltage. An important property of the res¬ 
onant exchange in graphene is that unlike non-resonant 
RKKY theory in graphene the type of the interaction 
(ferromagnetic or antiferromagnetic) is the same for the 
same-sublattiice (AA) and different sublattice (AB) con¬ 
figurations. The interaction is ferromagnetic on small 
distance if the bound state energy level lies in the con¬ 
duction band and antiferromagnetic for the bound state 
energy energy level lying in the valence band. Because 
the contributions from AA and AB configurations do not 
compensate each other in the resonant case, the random 
distributions of the magnetic adatoms results in a non¬ 
zero magnetic interaction. This finding opens a good 
possibility for studying the phenomena experimentally. 
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Appendix A: Calculation of eigenvalues 


Putting aside some complications introduced by 
graphene, our problem is the one of interaction be¬ 
tween two discrete levels and continuum, so-called Fano- 
Anderson problem 1 ^. The ’continuum’ can be repre¬ 
sented by a dense set of discrete levels, the determinant 
of the appropriate eigenvalue problem can be reduced to 
a 2x2 form: 
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The summation in ( A2|A3 ) is over all basis states, m,n = 
0,1, 2,.... The quantized momentum p in the tunnelling 
matrix element is discrete and corresponds to the index m 
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through the quantization condition |7]) . Indexes c, d cor¬ 
respond^© the electrons and holes, respectively. Equa¬ 


tion (A21 describes a shift of the adatoms energy levels 


due to interaction with graphene carriers. In the follow¬ 
ing calculation we change definition of the bound states 
energy levels as £, + A; -A £j. Note that the shift Ai 
is of the order of T 2 and can be neglected in our case 


. Equation (A3) describes interaction between adatoms 
via graphene electrons. Expanding (|Al|) we get: 


B = 


{e 1 +\ 1 -E){e 2 +\ 2 -E)-BB* =0, (A4) 
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where R = R 2 — Ri, k = E/hvp , fl u . c . = 3v3o 2 /2 is 
the graphene unit cell area, a is the C-C bond length, 
J ,, N t - are Bessel and Neumann functions of z-th order, 
respectively. The unit cell area emerges from ([ 5 ]) and 
(fl0| as L/N = Q u . c ./2nL. The roots of the characteris¬ 


tic polynomial (A4) can be parameterized with indexes 


E m ,i, where m describes the absolute value of momentum 
in unpe rturbed state E^, *=1,2 indicates two solutions 
of (A4) in accordance with the number of adatoms con¬ 


sidered (two). 


Appendix B: Exact equation for resonant indirect 
exchange interaction 


In this section we present an exact calculation of res¬ 


onant indirect exchange interaction following (A4). We 
obtained: 
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Let us consider the non-resonant limiting case. We as¬ 
sume (£0 = £1 = £ 2 ), £0 Ep, for undoped graphene 
we take Ep = 0. We analyse the two configurations: 
(i) adatoms are located on same sublattice Taa = t\A = 
t 2 A, tiB = t 2 B = 0, (ii) the adatoms are located on differ¬ 
ent graphene sublattice Tab = tiA = t 2 p, t 2A = tis = 0. 
For the two sublattice configurations we have 
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These results ( B3|B4 | describe indirect exchange inter¬ 
action in non-resonant limit and are in agreement with 
those obtained using RKKY theory in graphen d 14 * 15 iTEsI, 
the difference in the pre-factor resembles the specifics of 
our model where the bound state level exist also in a 
non-resonant case. Note that the absence of n factor in 
exact equation for g(E,r) (B2) is connected to the ir/2 
angle difference in our coordinate system froirP^ 17 . As 
one could see from (Bl) the indirect exchange oscillates 
with atomic length period. We make an averaging over 
these short-range oscillations: 


Z7T 

{E ex ) = j ^E ex { 7 ), 

0 

7= (K-K')R. 


(B5) 


This integral cannot be calculated analytically, but we 
can neglect the oscillation terms (i.e. cos( 7 ) -A 0) in 
(B2). We checked the error associated wi th d isregarding 
these terms by a numerical calculation of (B5) and found 
good agreement, the maximum error is less then 20% 18 . 
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